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ABSTRACT 

We explore the possibility of measuring the mass accretion rate (MAR) of galaxy clusters from their mass 
profiles beyond the virial radius /Quo■ We derive the accretion rate from the mass of a spherical shell whose 
inner radius is 2/Qoo. whose thickness changes with redshift, and whose infall velocity is assumed to be equal 
to the mean infall velocity of the spherical shells of dark matter halos extracted from iV-body simulations. This 
approximation is rather crude in hierarchical clustering scenarios where both smooth accretion and aggregation 
of smaller dark matter halos contribute to the mass accretion of clusters. Nevertheless, in the redshift range 
z = [0,2], our prescription returns an average MAR within 20 — 40% of the average rate derived from the 
merger trees of dark matter halos extracted from A ? -body simulations. The MAR of galaxy clusters has been 
the topic of numerous detailed numerical and theoretical investigations, but so far it has remained inaccessible 
to measurements in the real universe. Since the measurement of the mass profile of clusters beyond their virial 
radius can be performed with the caustic technique applied to dense redshift surveys of the cluster outer regions, 
our result suggests that measuring the mean MAR of a sample of galaxy clusters is actually feasible. We thus 
provide a new potential observational test of the cosmological and structure formation models. 

Keywords: galaxies: clusters: general - methods: A’-body simulations - cosmology: theory 


1. INTRODUCTION 

In the current model of the formation of cosmic struc¬ 
ture, where dark matter halos form from the aggregation 
of smaller halos, the mass accretion of dark matter halos 
is a stochastic process whose average behavior can be pre- 
dict ed with A’-body simulations and semi-analytical mod¬ 

els (Ivan den Boschll2002l:IZhao et al1l2003bt [Sheth & Tormen 
2004allbfc iGiocoli et alJl2007b IMcBride etakll2009l IZhao et al. 
200^ ). This process is generally investigated with the iden¬ 
tification of the merger trees of dark matter halos, en¬ 
abling the study of the mass accretion history (MAH) and 
the mass accretion ra t e (MAR) as a function of redshift z 
(Ivan den Boschl 120021 IFakhouri & Mai 120081 IMcBride et all 

l2009HFakhouri et al.ll2010t IGiocoli et al.ll2012l). 

Observationally, the exploration of the MAR of individ¬ 
ual dark matter halos has only been attempted on the scales 
of galaxies with the galaxy-galaxy merger rate: one usu¬ 
ally combines the number of observed pairs of close or dis¬ 
turbed gala xies with the theoretical merger probability and 
time scale dLotz et al.l 1201 It iJian et all 120 1 2b iCasteels et all 
120141) . Ho wever, different investigations reach discrepant 
conclusions dLotz et al .1120111) because the merger rate of dark 
mat ter halos is not ident i cal to the merger rate of galax¬ 
ies (IFakhouri & Mai 120081 iGuo & Whitel 120081 : iMoster et al.l 
1201 3l) and the two rates are related by dissipativ e processes 
which are difficult to model dHonkins et al.ll2013l) . 

In contrast, the measurement of the MAR of galaxy clusters 
can in principle be based on the estimated amount of mass 
in the cluster’s surrounding regions, where we can safely ne¬ 
glect the dissipative processes which affect the galaxy-galaxy 
merger rate. Nevertheless, no measurements of the MAR of 


galaxy clusters have been attempted so far. This observa¬ 
tional deficiency is due to the fact that in the large and less 
dense outer regions of clusters, galaxy members are diffi¬ 
cult to distinguish from foreground and background galaxies; 
other probes, e.g. X-ray emission, are below the sensitivity of 
current instruments. In addition, the outer regions of clusters 
are not in dynamical equilibrium and therefore the usual mass 
estimation methods based on virial equilibrium are inappro¬ 
priate. Considering this picture, one may think that the MAR 
predictions from iV-body simulations are not capable of being 
tested. 

Here we take a more optimistic perspective and explore the 
possibility of estimating the MAR of galaxy clusters by mea¬ 
suring the mass of a spherical shell surrounding the cluster. 
The thickness of this shell depends on the assumed infall time, 
on the radius at which the infall happens and on the initial ve¬ 
locity of the falling mass. This translates into a change of 
the shell thickness with redshift. Albeit rather crude when 
compared with the stochastic aggregation of dark matter halos 
in the hierarchical clustering formation model, this approach 
would provide a method to estimate the MAR that depends on 
the cluster mass profile at radii larger than the virial radius. 

In theoretical investigations, the relation between the mass 
density profiles of galaxy cl usters and their accretion history 
is known. For example, iLudlow et all (120131) find that the in¬ 
ner part of a halo retains the information on how the halo has 
accreted its mass through a correlation between the mean in¬ 
ner density within the scale radius r s and the critical density 
of the universe at the tim e when the mass of th e main progeni¬ 
tor is equal to M(< r s ). [Correa et al.l (I2015blfah confirm these 
findings and demonstrate that the MAH can be expressed with 
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a general formula _similar to the one originally proposed by 
iTasitsiom i et al.l (1 2004] ) for a ACDM model and widely stud¬ 
ied in iMcBride et al.l (120091) based on TV-body simulations: 
the MAH has an exponential evolution with redshift in the 
high 2 regime and follows a power law at low z when the 
accelerated expansion of the universe freezes the growth of 
pe rturbatio ns._ 

iDiemer & Kravtsov! (120141) find that the steepness of the 
slope of the outer halo density profile increases with in¬ 
creasing MAR. In addition, the central concentration is anti- 
correlated with the MAR. lAdhikari et al.l (120141) note that the 
location where the steepening of the slope is observed corre¬ 
sponds to the radius associated to the splashback of the ma¬ 
terial t hat the halo has recently accreted. Ivan den Bosch et all 
(120141) show that the growth of the central potential precedes 
the assembly of mass and introduce a formula that can be used 
to compute the average MAH in any ACDM cosmology with¬ 
out running numerical simulations. 

Here we derive a simple relation between the mass pro¬ 
file of a dark matter halo and its MAR derived from the halo 
merger tree extracted from TV-body simulations. This result 
is relevant because it implies that in principle we can esti¬ 
mate the MAR of galaxy clusters from the estimate of the 
mass profile in their outer regions. This mea surement can 
be performed with the caustic techniq ue ([Diaferio & Gelled 
1 1 9971 iDiaferioll 1999t ISerra et akll201 ll) which is not affected 
by the presence of substructures, the non-equilibrium state of 
the cluster outer region, an d the correlated large-scale struc¬ 
tures along the lin e of sight dDiaferioll 19991 ISeira et akll2Qlll 
iGeller etakll2013l) . The caustic technique only requires a suf¬ 
ficiently dense galaxy redshift survey in the cluster field of 
view to return a mass estimate accurate to 20% on average 
(ISerra et al.ll201 ll) . 

Estimating the MAR also requires the knowledge of the in¬ 
fall velocity of the shell. This quantity is inaccessible to ob¬ 
servations and remains a free parameter; we adopt the mean 
infall velocity of a shell surrounding dark matter halos in TV- 
body simulations as illustrated below. Although it prevents us 
from estimating the MAR of individual clusters, this choice 
allows us to cope with the unmeasurable velocity of the falling 
shell. Consequently, the technique we propose here aims at 
estimating the mean accretion rate of a sample of galaxy clus¬ 
ters rather than the accretion rate of individual clusters. 

We investigate the feasibility of our approach by con¬ 
sidering dark matter halos at redshift z < 2, with mass 
comparable to the clusters in the CIRS a nd HeCS catalogs 
dRines & Diaferioll2006l : iRines et mi2013 ~) whose outer mass 
profiles have already been measured with the caustic tech¬ 
nique. Our interest in investigating the growth of structures 
at nonlinear scales of galaxy clusters from an observational 
perspective, which is still relatively poorly explored ( e.g., 
iLemze et all 120131) . perfectly complements most of the cur¬ 
rent efforts that focus on constraining the growth factor in the 
linear and mildly nonlinear regimes with large-scale redshift 
surveys and weak-lensing tomography {e.g., Euclid, DES, 
eBOSS, DESI, PFS, LSST, and WFIRST). 

In Section [2] we introduce our spherical infall prescription. 
In Section [3] we discuss the prope rties of the CoDECS set of 
TV-body simulations (!Baldill2QI 2). We discuss our results in 
Section|4] 

2. THE ACCRETION RECIPE 

Our spherical infall prescription assumes a shell of mat¬ 
ter falling onto the enclosed halo. We aim to quantify M = 


dM/dt, where dM is the mass of a spherical shell of thick¬ 
ness Ar and proper radii Ri and Ri + Ar and dt is the time 
it takes to fall to Ri. Ri is the radius at which we con¬ 
sider the infall to happen. We choose A?’ = 5 s Ri, where 
S s is a free parameter. By solving the equation of motion 
d 2 r/df 2 = ao under the assumption of constant infall ac¬ 
celeration, ao = —GM(< Ri)/(Ri + 8 s Ri/2) 2 , and initial 
velocity Vi, we obtain the infall time ii n f from the equation 
aoff n f/2 + Viti n [ = —8 s Ri/ 2, where we consider the shell to 
be accreted when the shell middle point, initially at S s Ri/2, 
reaches Ri, we obtain 

t? ni GM{< Ri) - t in {2R 2 (l + S s /2fvi+ 

— RfS s (l + 5 a /2) 2 = 0 . 

We can now express the MAR obtained from our prescrip¬ 
tion as 

M = — = ^ 11 . (2 ) 

dt t in ( 

This recipe has three input parameters: the scale Ri that de¬ 
fines the infall radius, the initial velocity Vi and the thickness 
8 S of the shell. It is more convenient to use the infall time as 
an input parameter and to derive the shell thickness 8 S . In this 
case, Equation O reads: 

^ (^Ri + -^-Vitint^ + 

+ 8s (R$ + 2R 2 Vit in{ ) = 

GM(< Ri)tf n{ - 2 RjviUnt ■ 

Therefore, by choosing t ln t, we derive the thickness 6 s R, of 
the shell centered on the cluster and we can estimate the MAR 
of a galaxy cluster from its mass profile by using Equation (0. 

However, accretion is a stochastic process, and we need to 
verify that our simple approach is capable of correctly esti¬ 
mating the actual MAR. Below, we compare the MAR es¬ 
timated with our recipe with the MAR of dark matter halos 
derived from their halo merger trees obtained in an TV-body 
simulation. 

3. CODECS SIMULATIONS 

CoDECS (Coupled Dark Energy Cosmological Simula¬ 
tions) is a suite of TV-body simulations in different cosmo¬ 
logical models (see lBalclill20 1 2L for further details). Here we 
use the L-CoDECS simulation of a ACDM model, a colli¬ 
sionless TV-body simulation of a flat universe, with the follow¬ 
ing cosmological parameters consistent with the WMAP7 data 
(iKomatsu et al.l 12011 ): cosmological matter density 9 m o = 
0.226, cosmological constant Hao = 0.729, baryonic mass 
density ttbo = 0.0451, Hubble constant h = 0.703, power 
spectrum normalization erg = 0.809, and power spectrum in¬ 
dex n s = 0.966. The simulated box has a comoving volume 
of (1 Gpc ft -1 ) 3 containing (1024) 3 dark matter (CDM) par¬ 
ticles and the same amount of baryonic particles. The mass 
resolution is ?71dm = 5.84 x 10 10 M 0 h~ l for CDM parti¬ 
cles and mb = 1.17 x 10 10 M 0 /i _1 for baryonic particles. 
No hydrodynamics are included in the simulation. Baryonic 
particles are only included to account for the different forces 
acting on baryonic matter in the coupled quintessence models. 
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3.1. CoDECS Merger Trees 

Dark matter halos are identified at a given time with a 
Friends-of-Friends (FoF) algorithm with linking length b = 
0.2 times the mean interparticle separation. Flalos identi¬ 
fied with the FoF algorithm are called FoF halos. Flalo 
substructures are identified with the SUBFIND algorithm 
dSpringel et al.ll200l . We call SUBFIND halo the halo sub¬ 
structure that contributes the most particles to the FoF halo. 
To each FoF halo, we assign the radius fGoc0 and the corre¬ 
sponding enclosed mass M 200 of its SUBFIND halo. 

To derive the MAH of each halo, we use the 
merger trees provided in the CoDECS public database 
(www.marcobaldi.it/CoDECS). Each FoF halo at a 
given time G has a SUBFIND halo. We trace the main branch 
of this SUBFIND halo by searching for its main progenitor at 
each previous time step. The main progenitor of a SUBFIND 
halo at time ti is the SUBFIND halo at time t,_ 1 < G which 
contains the largest number of particles that will end up in the 
same FoF halo of the SUBFIND halo at time G- To each halo 
we associate the mass M200 of its SUBFIND halo and thus 
we derive the MAH as ( 1 / 200 ( 2 ) of the SUBFIND halos along 
the main branch. 

The definitions of the SUBFIND halo and the main progen¬ 
itor are not unique in the literature. The SUBFIND halo can 
also be defined as the substructure with the largest M 200 in the 
halo, rather than the substructure with the largest number of 
particles we adopt here. The usual definition of the main pro¬ 
genitor also is slightly different from ours: it can be defined 
either as the subhalo at G-i that donates the largest number 
of particles to the SUBFIND halo at G; or the most massive 
progenitor of the SUBFIND halo. In general, the differences 
between our MAHs and the MAHs obtained with the more 
common definitions of SUBFIND halo and main progenitor 
are negligible. However, our definitions guarantee that we al¬ 
ways trace the branch of the merger tree with the most massive 
halos. 

3.2. Cluster profiles 

We plan to apply our recipe for the estimate of the MAR 
to redshift surveys of galaxy cluster s by estimatin g their 
mass profile w i th the caustic met hod dDiaferio & Gellerl 1 997b 
lDiaferiolll999tlSerra et al.11201 11) . The caustic method returns 
mass profiles that are affected by 20-50% uncertainties when 
applied to clusters with M200 of 10 14 M 0 ft , -1 or larger. 
When applied to less massive clusters, the systematic errors 
introd uced by the caustic technique can become substantially 
larger dRines & Diaferioll20Tol; ISerra et al.ll20TTl) . Therefore, 
with this observational perspective in mind here we concen¬ 
trate on halos in the M 200 mass range 10 14 — 10 15 M 0 h~ l , 
which corresponds to the most common mass o f the clus¬ 
ters i n cata logs like CIRS dRines & Diaferiol2006l) and HeCS 
dRines et all 120131) . We thus define two mass bins at 2 = 0 
and follow the evolution with redshift of the clusters assigned 
to these bins. We select halos with median M 200 = 10 14 and 
10 15 M 0 h~ l at z = 0. The low-mass bin contains 2000 
objects at 2 = 0, while the high-mass bin is limited to 50 ob¬ 
jects: there are 36 clusters with M 200 > 10 15 M 0 h~ x in the 
simulated box, but we removed the 11 most massive clusters 
in order to have a mean mass similar to the median mass in 

9 f i/\ is the radius of the sphere of mass M& centered on a local minimum 
of the gravitational potential with average density A times the critical density 
of the universe, p c = SH 2 (z)/SnG, where H is the Hubble parameter. 


the bin. Table |T] lists the mean, standard deviation, median, 
68%, and 90% percentiles of each mass bin at z = 0. 

For each halo in the two mass bins an d fo r each progenitor 
at higher redshift (as defined in Section I3TT1 ). we evaluate the 
mass profile and the profile of the radial velocity u r;l <i up to 
10 f? 2 oo- We show the radial velocity profile in the two bins at 
different redshifts in Figures]]] and |2] 

With the radial velocity profile we can identify three re¬ 
gions: an internal region with v rm \ ~ 0, where matter is or¬ 
biting around the center of the cluster; an infall region, where 
u r ad becomes negative and indicates an actual infall of matter 
toward the center of the cluster; and a Hubble region at very 
large radii, where u ra( j becomes positive and the Hubble flow 
dominates. Broadly speaking, the infall radius R ln f, i.e., the 
radius where the minimum of u ra d occurs, is between 2/Goo 
and 3i?2oo. in depend ent of mass and redshift. 

[More et ahl (120151) suggest the use of the splashback radius 
as the physical halo boundary because it separates the infall 
region from the region where the matter has already been ac¬ 
creted. The splashback radius is defined as the outermost ra¬ 
dius reached by accreted material in its first orbit around the 
cluster center. Its exact location depends both on redshift and 
MAR, but it is in general larger than /Goo- Noticeably, they 
also show that the splashback radius is not affected by the evo¬ 
lution of the critical density of the universe, unlike the usual 
/Guo, whose definition is based on the average overdensity of 
the halo at a given redshift. By adopting /Goo, part of the evo¬ 
lution of the halo properties with redshift simply is a conse¬ 
quence of the evolution of the critical density of the universe; 
this effect gener ates a so-ca lled pseudo-evolution of the halo 
properties (see iDiemer et all 120131) . This pseudo-evolution 
substantially disappears when we adopt the splashback radius. 
This use of this radius might thus be more preferable in the in¬ 
vestigations of the redshift evolution of the halo properties. 

Our resultsju Fi gures [H and [2] are consistent with the re¬ 
sults of lMore et al.l (120151) 7 For massive objects at z — 0, the 
splashback radius R sp is close to 2 /Goo and the infall radius 
is i?i n f ~ 1.4i? sp . Given that M(< R sp ) and M(< i?inf) 
are not affected by the pseudo-evolution mentioned above, we 
use Ri = 2 f ?200 ~ Rsp as the radius at which we consider 
the infall to happen in our spherical infall prescription and 
Vi = Usheii ~ ^rad(-Rinf) (see Equation ®). 

For each halo at 2 = 0 in the two mass bins we build 
the MAH at 2/Goo (Figure |3}. In Table [2] we list the mean, 
standard deviation, median, 68%, and 90% percentiles of 
M(< 2f?2oo) at 2 = 0 for the tw o mass b ins. In Fig ure [3] 
we also show the MAH model by iGiocoli et all (120130 (see 
Appendix) rescaled to 2 /Goo- We obtain this rescaling by 
extending the mass density profile to 2/Goo using the NFW 
functional and adopting the redshift-dependent relation be¬ 
tween the concentration and mass of |Zhao et al.l (l2003al) mod¬ 
ified for M200 bv lGiocoli et all (12013b . For a direct compari¬ 
son with previous results in the literature we also estimate the 
MAH at the standard /Goo- We report this analysis and the 
comparison with previous work in the Appendix. 

We find the r atio M{< 2iGon) MGno to be ~ 1.6, similar 
to the result of lMore et alJ (120151) for M sp /M 200 in massive 
objects at 2 = 0. In real observations measuring the infall 
radius of a cluster, namely where the infall velocity reaches 
its minimum is currently unfeasible so we keep fixed R,, = 
2R 200 f°r all masses at all redshifts. This choice clearly is 
an oversimplification, but Figures [T] and [2] indicate that this 
assumption is reasonable. Indeed J?j n f lies within 2/Goo and 
3fGoo for a wide range of masses and redshifts. 
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[10 14 M 0 h - ’] (2000) 



[10 14 M 0 h - '] (2000) 



r /^200 


Figure 1. Radial velocity profile for clusters in the 10 14 Mg h ~ 1 mass bin at z = 0 (left panel) and 2 = 1 (right panel). The thick-solid and thick-dashed 
lines indicate the mean and median profiles, respectively. The thin-solid and thin-dashed lines indicate the standard deviation and the 68% of the distribution, 
respectively. 


[10 15 M© h" 1 ] (50) 



f /^200 


[10 15 M© h" 1 ] (50) 



r /^200 


Figure 2. Same as Figure^except for the 10 15 Mq h 1 mass bin. 

[10 14 M© h' 1 ] (2000) [10 15 M© h' 1 ] (50) 




Figure 3. MAH at 2 R 200 for the 10 14 Mq h ~ 1 (left panel) and 10 15 Mq h ~ 1 (right panel) mass bins. The thick-solid and thick-dashed lines indicate the 
mean and median MAH, respectively. The thin-solid and thin-dashed lines indicate the standard deviation and the 68% of the distribution, respectively. The 
t hin-dotted-dashed line shows indicates 90% of the distribution. The blue curve indicates the mean MAH for M 2 oo- The green curve indicates the model by 
IGiocoli et alj 12013b rescaled to 2 R 200 . 
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Table 1 

Mean, Median, Standard Deviation, and Percentile Ranges of M 200 and Number of Halos in the Two Mass Bins at 2 : = 0 


Mean 

a 

Median 


68% 


90% 

Number of Halos 



M 200 [Mq h 

- 1 ] 




1.00 X 10 14 

2.90 x 10 12 

1.00 x 10 14 

(0.97 

- 1.04) x 10 14 

(0.96 - 

1.05) x 10 14 

2000 

1.04 X 10 15 

1.26 x 10 14 

1.00 x 10 15 

(0.91 

- 1.19) x 10 15 

(0.88 - 

1.24) x 10 15 

50 


Table 2 

Mean, Median, Standard Deviation, and Percentile Ranges of M(< 27^200) in Each Mass Bin at z = 0 


Mean 

C7 

Median 

68% 

90% 

M(< 2 R 200 ) [M 0 h- 1 ] 

1.60 X 10 14 

1.73 x 10 13 

1.58 x 10 14 

(1.45 - 1.74) x 10 14 

(1.37- 1.90) x 10 14 

1.55 X 10 15 

2.41 x 10 14 

1.56 x 10 15 

(1.32 - 1.82) x 10 15 

(1.19 - 2.00) x 10 15 
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4. OUR ACCRETION RECIPE VERSUS THE MERGER TREES 

When we observe a cluster we look at a particular instant of 
its evolution and we can trace the evolution of its mass neither 
backward nor forward in time. Therefore, the MAH is not a 
measurable quantity. In contrast, the MAR can in principle 
be estimated as the ratio of the mass that is being accreted 
by the cluster at a given time and its infall time as described 
by our spherical infall prescription in Section[2] Equation 0 . 
This estimate depends on the external mass profile and on the 
mean infall velocity. The mass profile is measurable with the 
caustic method, whereas the initial velocity i\ (see Equation 
( 0 ) remains a free parameter which can be inferred from N- 
body simulations. 

As already stated in Section 1-01 we choose Ri = 2f?2oo ~ 
i? sp as the radius where the infall occurs. As initial velocity v, 
for the infall, for each redshift in each mass bin we adopt the 
mean velocity in our first radial bin [2 — 2.5]f?200 ~ R'mf 
(see Figures Q] and 0: this velocity ranges from —200 to 
—250 km s -1 . We call this velocity u s heii- The value of u ra d 
in the radial bin [2.5 — 3]f?2oo roughly remains in the same 
range and adopting this radial bin instead of the first bin has 
negligible effects on the final results. Although our prescrip¬ 
tion for the MAR estimate clearly depends on the choice of 
both Ri and u s h e ii, using the mean infall velocity instead of 
the infall velocity of each halo enables the design of a feasi¬ 
ble procedure for the observational estimate of the MAR of 
clusters. 

The spread of the distribution of v s h e ii, estimated as a = 
a Vi /VN, where N is the number of halos used to evalu¬ 
ate '(.’shell, is smaller than the spread of the distribution of 
the individual Vi s, a Vi , shown in Figures Q] and [2] For the 
10 14 Mq/i _1 mass bin, the la relative spread of u s heii is 
smaller than 2% and propagates into a relative spread on the 
MAR of 5% at most, which is well within the 68% percentiles 
of the MAR distribution obtained with u s heii- Similarly, for 
the 1O 15 M 0 h 1 mass bin, where N drops from 2000 to 50, 
the lcr relative spread of u s heii is ~ 14%, implying a ~ 35% 
relative spread on the MAR. Incidentally, even by using the 
lcr spread of the distribution of the individual Vi s, the result¬ 
ing MAR’s are within the 68% percentiles of the distribution 
of the MAR’s obtained with u s heii for both mass bins. The dif¬ 
ferent spread of v s heii in the two mass shells deriving from the 
different size of the two cluster samples indicates that the ex¬ 
pected uncertainty on the MAR estimated with real data will 
vary substantially depending on the number of observed clus¬ 
ters, as one can expect. We plan to fully quantify this un¬ 
certainty in future work by estimating the MAR of simulated 
clusters from mock redshift surveys. 

For the infall time, which is the last parameter of the model, 
we choose t ln f = 10 9 yr. This value is suggested by the 
redshift-independent ~ 1 Gyr time step of the snapshots of 
the simulation used to estimate the MAR from the merger 
trees; it also has the advantage of being similar to the dynami¬ 
cal time, simply defined as f<jyn ~ R/cr, for the clusters of our 
analysis. In fact, for the 10 14 —10 15 M 0 h~ l clusters at z = 0 
we consider here, R ~ 1 Mpc and cr ~ 1000 km s _1 , and 
fdyn ~ 10 9 yr. For the progenitors of these clusters at higher 
redshifts, which have masses at most a factor 10 smaller, the 
velocity dispersion a is If) 1 / 2 ~ 3 times smaller and the virial 
radius R is smaller by roughly a similar factor. Therefore, 
fdyn remains basically constant and equal to 10 9 yr. Finally, 
this equality between t<jyn and the snapshot time interval pre¬ 
vents us from assuming t - ln f very different from 1 Gyr; if t- ln { 


departs too much from the snapshot time interval, the compar¬ 
ison of our estimated MAR’s with the MAR’s extracted from 
the merger trees would be inappropriate. We therefore inves¬ 
tigate the dependence of our results on t\ n f within ~ 20% of 
10 9 yr and find that our results remain unaffected. 

Once Ri, u s heii and t; n f are specified, the model is com¬ 
pletely determined by Equation 0 . For each halo in the two 
mass bins and for each progenitor at higher redshift, we eval¬ 
uate the thickness S s of the infalling shell and its mass. We 
show the evolution with redshift of the shell thickness S s in 
Figure ® The shell thickness increases with increasing red¬ 
shift, and the intrinsic scatter of the distribution also enlarges. 
This fact reflects the individual evolution of M(< 2f?2oo), as 
shown in Figure [3] The solid blue line in Figure [4] marks the 
value S s = 0.5 for which the external radius of the shell is 
equal to 3f?200, close to the cluster turnaround radius. This 
value of S s is reached between redshift 1 and 1.5, depending 
on the mass of the cluster. 

The accretion onto a cluster is a highly anisotropic process; 
nevertheless, we are confident that, given the thickness of the 
shell, we are taking into account almost all of the mass that is 
actually falling in the time interval £j n f. 

Figure [5] shows the MAR of all the clusters in the two mass 
bins estimated with Equation 0 . It also shows the MAR de¬ 
rived from the merger trees of the halos. 

We see that the mean and median results from the merger 
trees lie within the region defined by the 68% percentile range 
of the distribution of the MAR’s obtained from Equation 0 
for both mass bins. For the 10 15 M 0 /i _1 bin, the mean and 
median values of our estimate are ~ 40% smaller than the 
average MAR values from the merger trees. This systematic 
underestimate might be due to a statistical fluctuation because 
the sample only contains 50 clusters compared with the 2000 
clusters of the less massive bin. In fact, a similar underesti¬ 
mate is observed when we compare the MAH from the merger 
trees with the Giocoli model in the right panel of Figure 0 In 
contrast, for the 10 14 M 0 /i _1 mass bin, the mean and me¬ 
dian MAR from our prescription recover the merger tree re¬ 
sults within 20% in the redshift range 2 = [0,2]. Our results 
are relevant because they show that our simple spherical in¬ 
fall prescription can in principle provide a method to estimate 
the average MAR of galaxy clusters from redshift surveys. In 
future work, we will apply our prescription to synthetic red¬ 
shift surveys of galaxy clusters to quantify the uncertainties 
and possible systematic errors of our procedure. 

Clearly Figure 0only compares the average MAR obtained 
from the merger trees of individual halos with the average 
MAR provided by our spherical infall technique. Our recipe 
was not conceived to completely capture all the features of the 
MAR derived by the complex merging process of individual 
halos. Nevertheless, the average of the MAR of individual 
halos still is satisfactorily estimated by our recipe. 

Figure [6] shows the distribution of the ratio between the 
MAR estimated with our recipe and the MAR derived from 
the merger tree for each individual halo of the 10 14 M 0 /i _1 
mass bin at four different redshifts. The distribution has a tail 
toward large values. Remarkably the median value of this ra¬ 
tio is close to the ratio between the average MAR from our 
model and the MAR from the merger trees (red line). The 
mean clearly is larger because of the tail of high values. The 
68% percentile ranges from 0.2 to 2.2. This result confirms 
that with our approach we can measure the mean MAR but 
not the MAR of the individual clusters, which is affected by 
lack of spherical symmetry and large variations of the infall 
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Figure 4. Redshift evolution of the shell thickness & s for the 10 14 Mg h ~ 1 (left panel) and 10 15 Mg h ” 1 (right panel). The thick-solid and thick-dashed lines 
indicate the mean and median <5 S , respectively. The thin-solid and thin-dashed lines indicate the standard deviation and the 68% of the distribution, respectively. 
The solid blue line marks the value S s = 0.5 for which the external radius of the shell is equal to 3f?200- 




1 2 3 1 2 3 

1+z 1+z 


Figure 5. Results of our spherical infall model for the 10 14 Mg h~ 1 (left panel) and 10 15 Mg (r 1 (right panel) mass bins and comparison with the MAR 
from merger trees. The blue solid and dashed lines are the mean and median MAR from Equation j2}. The green area indicates the 68% of this MAR distribution. 
The mean and median MAR from the merger trees are indicated by the black solid and dashed lines. Residuals from the median and the mean values are shown 
in the insets at the bottom of each panel. 


velocity. 

Figures [T] [2] and 0 show the intrinsic halo-by-halo scatter 
both in radial velocity and mass. Our choice to use the same 
t> s heii for all clusters in a given mass bin at a given redshift 
implies that we neglect the scatter that originates the spread 
of the distribution shown in Figure [6] However, this choice 
keeps the model relatively simple and applicable to real clus¬ 
ters. It is worth saying that the impact of the large-value tail 
is reduced if we take the ratio of the averages of the MAR of 
each individual halo estimated with our prescription and with 
the merger trees, rather than the average of the ratio. The ra¬ 
tio of averages is shown with the red lines in Figure [6] and it 
corresponds to the result shown in Figure 0 The remarkable 
and encouraging result of our analysis is that the agreement 
shown in Figures [5] and [6] is obtained without requiring any 
input information from the merger trees. 

5. CONCLUSIONS 

We investigate the feasibility of directly measuring the 
mean MAR of a sample of galaxy clusters from their mass 
profile. To measure the mean MAR we suggest a prescription 


based on the infall of a spherical shell of matter onto the halo, 
with constant acceleration and initial velocity derived from 
the average infall velocity of matter around dark matter halos 
in 2V-body simulations. Once we fix the scale Ri, that defines 
the halo radius at which the infall occurs, the initial velocity 
Wsheii and the infall time ij n f our method only depends on the 
mass profile at large radii, beyond Ri. 

We consider dark matter halos from the CoDECS set of N- 
body simulations and compare their MAR estimated with our 
prescription with the MAR estimated from the merger trees 
extracted from the simulations. We focus on two sets of halos 
with mass M 2 oo around 10 14 and 10 15 M 0 h~ r at z = 0. 

We recover the mean and median MAR obtained from the 
merger trees without bias and with 20% accuracy in the red¬ 
shift range z = [0, 2] for the 10 14 M 0 h~ x mass bin. The 
accuracy is about 40%, with a systematic underestimation of 
~ 40%, for the 10 15 M 0 h~ x mass bin. This result is im¬ 
pressive given the simple assumptions of our prescription and 
the fact that no input parameter of the model is taken from the 
merger trees. Our result does show that measuring the mean 
MAR of a sample of real galaxy clusters is in principle possi- 
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Figure 6. Histogram of the halo-by-halo ratio between the MAR from our spherical infall model and the MAR from the merger trees at z = 0 (top-left panel), 
z = 0.35 (top-right panel), z = 0.55 (bottom-left panel) and 2 = 2 (bottom-right panel) for the 10 14 Mg h ~ 1 mass bin. The solid, dashed and dotted-dashed 
vertical lines mark the mean, median and 68% of the distribution, respectively. The red line marks the ratio between the average MAR from our model and the 
average MAR form the merger trees. 


ble. 

A fundamental step to assess the feasibility of our approach 
is to apply the caustic technique to realistic mock redshift sur¬ 
veys of galaxy clusters extracted from TV-body simulations 
and quantify the accuracy of the estimated MAR. This inves¬ 
tigation remains for further studies. Similarly the analysis of 
the dependence of the method parameters and its results on 
the cosmological model and on the theory of gravity remain 
to be investigated. Specifically, i> s heii might turn out to vary 
substantially with the assumed cosmology and thus to be a 
crucial parameter of the method. 

We might expect that the accuracy will be better than 20% 
on the average MAR if we estimate the MAR of individual 
clusters in a given mass bin and estimate their average MAR; 
this approach would minimize the systematic errors due to 
projection effects which dominate th e estima te of the mass 
profile with the caustic technique. [Rines et all (2013) have al¬ 
ready applied this approach to measure the total mass of clus¬ 
ters within their turnaround radius, the ultimate mass M ta . By 
combining 50 CIRS clusters with 58 HeCS clusters, they find 
M ta /M 2 oo = 1-99 ± 0.11, a measure accurate to 5% and in 
agreement with the ACDM prediction where M ta /M 2 oo has 
a log-normal distr ibuti on with apeak at mass ratio 2.2 and 
dispersion 0.38 (iBusha et al.ll2005i) . 

Our measurements of the average MAR may provide an 
additional tool to discriminate among different cosmological 
models if deviations from the ACDM model generate mea¬ 


surable differences in the MAR. 
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APPENDIX 

In this appendix, we investigate the MAH of M 200 of dark 
matter halos: we compare the merger tree results with two 
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known fitting formulae and a theoretical model. We calculate 
the mean and median MAH for the objects in four mass bins 
(Table [3} The two most massive bins are the same used in 
the main body of the paper. The four bins are centered on 
M200 = 10 15 ,10 14 ,10 13 , 5 x 10 12 M 0 ft- 1 . Each bin con¬ 
tains 2000 objects at z = 0 with the exception of the largest- 
mass bin whose sample is limited to 50 objects. The lowest- 
mass bin is centered on M 200 = 5 x 10 12 M 0 ft -1 , because 
M 200 = 10 12 M 0 ft -1 is below the resolution limit of about 
100 particles per subhalo set by SUBFIND. Table[3]also lists 
the 90% percentile mass range of each mass bin. 

Figure [7] shows the MAHs of the four mass bins. We limit 
our study to the low redshift range z = [ 0 , 2 ] because we are 
interested in the observational relevance of our analysis. As 
we can see in Figure[7] for the two largest-mass bins the mean 
MAH agrees with the median MAH within 20%. In the two 
smallest-mass bins the difference between the mean and the 
median MAH is never larger than 5%. In all four cases the 
standard deviation and the 68 % percentiles are comparable. 
The results from the largest-mass bin are noisier because of 
the low-number statistics. The number of objects A'f ,1 at each 
z decreases with increasing z, due to the resolution limit: not 
all the objects selected at z = 0 have merger trees that reach 
z = 2. Indeed, the decrease is larger for less massive objects 
which are already closer to the mass resolution limit at z = 0 . 

A.l. Fitting Formulae 

Different fitting formulae for the MAH shown in Figure [7] 
exist in the literature. We focus on two of them. B y usi ng the 
extendedPress-Schechter (EPS) formalism dBond et a l.ll 1991 h 
lEacev & Colell 1 99~3l) . lvan den Boschl (120021) proposed 


M(z) = M 0 exp 



ln(l + z) 

ln(l + z/) 



(4) 


where Zf and v are free parameters. The redshift formation 
Zf indicates the redshift when the halo has accreted half of 
its final mass. We fit both the median and the mean MAH 
with the equation above by assuming Poisson errors weighted 
by the number of halos in each redshift bin. We quantify the 
deviation from this analytic description of the MAH with the 
rms of the fit 


(r.m.s.) 2 


1 | _M(zi) — M mo( jel(Zi)] 


(5) 


We list the best-fit parameters of Equation (0i along with 
the rms of the fits in Table [4] As expected in hierarchical 
clustering scenarios, the value of the best-fit parameter Zf in¬ 
creases with decreasing mass because more massive objects 
tend to form later than less massive ones. 

The second formula we considered was fir st proposed by 
iTasitsiomi et alJ (120041) and widely studied by McBri de et al.1 
(120091) : 


M(z) = M 0 (l + zfe~^ z , ( 6 ) 

where ft and 7 (> 0) are free parameters. This formula repre¬ 
sents an exponential growth with a redshift-dependent correc¬ 
tion. It is a revision of the sim ple one-pa rame ter exponential 
form M(z) = Moe~ az (iWechsler et all 120021) . where a = 
ln( 2 )/Zf. By using the EPS formalism. ICorrea et akl(l2015al) 
showed that in a ACDM model the exponential growth is a 
good description of the MAH at high z, while the power-law 


behavior at low z is necessary because the accelerated expan¬ 
sion of the universe slows down the accretion. For this reason 
Equation (( 6 )) appears to be a general description of the MAH 
of dark matter halos in a ACDM model independently of the 
halo mass. The param eters ft and 7 are related to the power 
spectrum (ICorrea et al.ll2015al) . The value of ft — 7 is a param¬ 
eter describing the mass growth rate at low redshift. We use 
Equation 0 to perform a fit with Poisson errors weighted by 
the number of halos in each redshift bin and evaluate the rms 
as in Equation (0. We list the best-fit parameters of Equation 
([0) and the rms of the fits in Table [4] 


A.2. Comparison with a Theoretical Model 

By follow ing and generali zing the form ali sm o f 
Lacev&Cold (fl993l) and iNusser & Shetfl (fl999h . 
Giocoli et al.l (12012b introduced a new theoretical model 
to describe the MAH of dark matter halos. This simple 
model, which enables the derivation of a generalized redshift 
formation distribution, has already been applied to the 
CoDECS simulations in lGiocoli et all (120131) . 

Here we su mmar ize the relevant definitions and refer to 
iGiocoli et al.l (120121) for further details. The model defines 
the redshift formation Zf of a halo of a given mass Mq at a 
given redshift zq as the redshift when the object has accreted 
a fraction / of its final mass Mq, for any fraction o</< 1 . 
The variance of the linear fluctuation field is 


S(M) = A 

where W(kR) is a top-hat window function of scale R = 
(SM/Anpm) 1 / 3 , Pm is the comoving background density of 
the universe, and Pn n (k) is the linear power spectrum. 

The initial threshold overdensity for spherical collapse is 


W 2 {kR)P\m{k)k 2 dk 


(7) 


Sc{z) 


^c,DM 

D^z)’ 


( 8 ) 


where <5 c ,dm is the linear overdensity at redshift z and D+(z) 
is the growth factor normalized to unity at z = 0. 

The cumulative formation redshift distribution, in terms of 
the scaled variable 


is given by 


S c (zf) - S c {zp) 

\/S(fM) - S(M) ’ 


(9) 


p {> w f) 


a f 


w 2 f /2 . -1 

e f/ + a/ — 1 


( 10 ) 


The model has a single free parameter a / which depends 
on the fraction / of the final mass that is used to define the 
formation redshift z t. Fo r the same set of simulations we use 
here. IGiocoli et al.l ( 2013 ) find 


1.365 _ 2 f3 

a f = Jo*5 e 


( 11 ) 


Since Equation ( ITOl) can be inverted, it is possible to eval¬ 
uate the median redshift Zf when a halo accretes a fraction / 
of its final mass M with the relation 


Sc(zf) = $c(zq) + WfVStfM) - S(M) , ( 12 ) 
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Table 3 

Mean, Median, Standard Deviation, and Percentile Ranges of M 200 and Number of Halos in Each Mass Bin at z = 0 


Mean 

a 

Median 


68% 


90% 

Number of Halos 



M 200 [Mq h- 

- 1 ] 




1.04 x 10 15 

1.26 x 10 14 

1.00 x 10 15 

(0.91 

- 1.19) x 10 15 

(0.88 - 

1.24) x 10 15 

50 

1.00 x 10 14 

2.90 x 10 12 

1.00 x 10 14 

(0.97 

- 1.04) x 10 14 

(0.96 - 

1.05) x 10 14 

2000 

1.00 x 10 13 

2.36 x 10 10 

1.00 x 10 13 

(1.00 

- 1.00) x 10 13 

(1.00 - 

1.00) x 10 13 

2000 

5.00 x 10 12 

6.04 x 10 9 

5.01 x 10 12 

(5.00 

- 5.01) x 10 12 

(5.00 - 

5.01) x 10 12 

2000 


Table 4 

Best-fit Parameters of the Formulae by Van Den Bosch, Equation 0, and McBride, Equation 0 




Median 




Mean 







van den Bosch 




M 2O o[1O 14 M 0 h- 1 } 

z f 

V 


rms 

z f 

V 


rms 

10 

0.358 ±0.011 

1.227 ±0.035 


0.044 

0.381 ± 0.012 

1.252 ±0.037 


0.043 

1 

0.711 ±0.003 

1.702 ±0.012 


0.006 

0.709 ± 0.003 

1.582 ±0.012 


0.022 

0.1 

0.944 ± 0.005 

1.675 ±0.019 


0.025 

0.925 ± 0.005 

1.519 ±0.018 


0.015 

0.05 

1.033 ±0.010 

1.317 ±0.024 


0.032 

1.036 ±0.011 

1.242 ±0.022 


0.024 





McBride 




M 2 OO [io 14 M 0 h- 1 ] 

P 

7 

P ~1 

rms 

P 

7 

P -1 

rms 

10 

-1.000 ±0.270 

1.150 ±0.180 

-2.160 

0.048 

-0.690 ±0.270 

1.280 ±0.180 

-1.970 

0.028 

1 

±1.283 ±0.044 

1.920 ±0.029 

-0.637 

0.030 

±0.812 ±0.043 

1.566 ±0.028 

-0.754 

0.010 

0.1 

±0.898 ±0.051 

1.342 ±0.035 

-0.444 

0.039 

±0.480 ±0.051 

1.065 ±0.035 

-0.585 

0.029 

0.05 

-0.005 ± 0.070 

0.658 ±0.051 

-0.664 

0.053 

-0.217 ±0.069 

0.513 ±0.050 

-0.730 

0.037 


where 

Wf = yj 2 ln(a/ + 1) (13) 

is the median value of vu defined by the usual relation P(> 
u>f) = 1/2. Equation (fT2l> can be translated into a MAH for 
a given final mass M 0 . 

We compare the simulation results and the Giocoli model in 
Figure [7] The model of iGiocoli et al.1 (120 1 2l) is built by eval¬ 
uating the median redshift at which the halo has accreted a 
fixed fraction of its final mass whereas we evaluate the mean 
and the median M 200 for all the objects in a given mass bin 
at each redshift. We list the rms defined in Equation <|5]» in 
Table [5] In general, the global agree ment between the the¬ 
oretical model of IGiocoli et al.l (120121) and the simulation in 
each mass bin is similar for both the mean and the median, as 
it can be seen from the rms values in Table 0 In the largest- 


mass bin the model overestimates the MAH obtained from 
the simulation, while in the other three bins the agreement is 
within a few percent up to z ~ 1. Toward higher z, the model 
starts underestimating the simulation MAH. This discrepancy 
is more pronounced and appears at decreasing redshifts for 
decreasing halo mass. This behavior originates from the mass 
resolution and the consequent decrease of the number of halos 
iVh a i at a given redshift. 
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